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ABSTRACT 


The successful application of holographic interferometry, and an 
associated mathematical reduction process , to the determination of an 
asymmetric three-dimensional density field of an aerodynamic phenom- 
enon is reported. 

An integra] inversion method from the field of plasma physics has 
been computerized, extensively evaluated and applied to the determin- 
ation of functions , both axisymmetric and asymmetric, which simulate 
aerodynamic density fields. 

The application of holographic interferometry has been extended 
to provide multiple holograms about a test region, with sufficient 
coverage to provide interferometric data for the successful solution 
of the density field. 

The analytical and experimental methods developed were applied 
to an experimental axisymmetric test field, the supersonic flow from 
a free jet, and shown to be comparable to a previous solution obtained 
by the Abel inversion method. Further, the free jet was tilted to 
provide a test field which was asymmetric in the plane of solution. 
Comparison of the resulting asymmetric solution was shown to be self- 


consistent with the previously obtained axisymmetric solution. 
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I. INTRODUCTION 


Holography has enabled the aerodynamicist to "freeze" the interfero- 
metric view of a transient phenomenon and subsequently to view the field 
in three dimensions, as through a "window." This new capability has 
indicated the possibility of quantitatively determining the density in such 
a field, with no symmetric restrictions on its form. 

To make quantitative determinations of asymmetric density fields, 
it was necessary to: 


1) Invert the fringe number functions that describe the interfero- 


x. 


3 9) = ¥ f(x,y) dx’ 


to obtain the asymmetric density function f from within the integral. 


metric data, 


In the search for such an inversion process the possibility of using an 
asymmetric line integral inversion scheme which had been developed in the 
field of plasma emissivity [Maldonado, 1966], became evident. The method 
was evaluated by extensive computerized testing on density functions 
found in aerodynamics, including shock waves. The inversion method was 
found to be quite accurate, giving results which were generally within one 
percent of the test function. 

2) Obtain holographic interferograms of sufficient angular coverage 
about a test region to provide the interferometric data required by the in- 


version process. To obtain the necessary experimental data, a holographic 





work table was pore cricted and used to design the necessary optical 
arrangement. The capability of achieving multiple Q-switched holographic 
interferograms about a test region, with no opaque objects present, was 
demonstrated. 

3) Apply real holographic Per rercictric data to an inversion 
process in order to demonstrate the practicality of the technique. A 
stepwise evaluation of the system included an axisymmetric evaluation 
of a free jet for comparison with the previously available solution of 
Winckler [1948]. Having attained a solution for the three-dimensional 
axisymmetric field, the asymmetric capability of the technique was 
tested by tilting the free jet, which destroyed the axisymmetry in the 
plane of inversion. The resulting data were inverted asymmetrically 
and compared with the axisymmetrically measured density. The result 


is shown to be a self-consistent comparison of the asymmetric inversion. 


II. THEORETRICAL ANALYSIS 


A. THE INTERFEROMETRY OF AERODYNAMIC FIELDS 

Holography has opened a new branch of interferometry for aero- 
dynamic research. Ordinary Mach-Zehnder type interferograms have 
been restricted to the analysis of two-dimensional or axisymmetric flow 
fields. Such interferograms are able to examine only the light beam which 
has passed through the test region from a single direction. Holograms are 


able to reproduce a light wave in both amplitude and phase, and further to 





simultaneously reproduce two such light waves which have occurred at 
different times. The resulting comparison of two reproduced light waves 
yields an interferogram which allows the analysis of the light which has 
passed through a test region from many different directions. 

By interfering two coherent light waves with each other, interfero- 
metry compares their phases, and thus provides a measured comparison 
of their integrated index of refraction over their respective distances 
of travel from a common source to the film. This relative phase, 
expressed in multiples of wavelengths and commonly called fringe num- 


ber, is described by the equation: 


« 


Le 
9 = < \ [n (x,Y,2) -n,(x,v,2) dh (1) 


where: X = the wavelength of the light at the film, 
n =the index of refraction seen by the traveling light wave, 
n, = the index of refraction seen by the wave of the comparison 
scene , 


L-L = the total snes of travel for each of the waves. 

The index of refraction of a gas is well represented as a function of 
the gas density [Liepmann and Roshko, 1957] by the first two terms of a 
Taylor expansion: 

Beta SF xs (2} 
where: Ne = the gas density at which @ is measured, 


@ =the measured coefficient. 


For air at 0°C. and 760 mm. He t> = ,0012929 gm/cc. , and for deep 
red light ,@= .000291. Thus one can rewrite equation (1): 
g (* 
— 93> f (K,y,2) -& (,v,2)| dX (3) 
XK N Ss 
Equation (3) is a line integral equation for the unknown density. With 
sufficient evaluation of the function g (a function of three spatial 
variables which one obtains from the interferograms), one can invert the 
equation to obtain the density function. 
1, The Mach-Zehnder Interferometer 
Figure (1) is the schematic arrangement of a typical Mach-Zehnder 
interferometer [Ladenberg, ed. 1954]. The two light waves which are 
compared travel paths A and B respectively. Alignment of the instrument 
is critical to within the order of a wavelength of light, and the mirrors 
and lenses must be distortion free to one-tenth of that value. The chief 


practical disadvantage , however, is that one can evaluate the light beam 
from only one direction at a time through the test section. 

By properly aligning the components with a uniform density field 
in the eee section, one can exactly recombine the two light beams A and B 
at the film plane of the camera. Then with a flow in the test section, 
any phase differences recorded are those due to variations in the test 
section density. Such an interferogram is called "infinite fringe” 

If one purposely misaligns one beam with the other (by rotating a 


mirror slightly) a deliberate fringe pattern is recorded. Variation in test 


Section density is then recorded as a variation in the regularity of the 
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fringe pattern. This second variety of interferogram is called "finite 
fringe." 
2. Holography 

First reported by Denis Gabor [ 1948, 1949, 1951 ], holography is 
itself an interferometric process, By combining two coherent beams 
directly on a photographic film (Figures 2, 2a), one records a micro- 
scopically small interference pattern. Positions on the film where the 
two waves are in phase will record darkly. Positions where they are out 
of phase will remain light. When the developed photographic plate is 
reilluminated by the reference beam (Figures 3, 3a), the interference 
Me cen acts as a very complex diffraction grating and FEE eEACee three 
beams; these are the zero order, or transmitted illumination beam, and 
the positive and negative first orders, which are re-creations of the 
original object beam. 

Mathematically, one ean represent the phenomenon as the com- 
bination of any two general waves which obey the scalar wave equation 


[Brandt, 1968]. Let a wave be described as follows: 


iE At) = A(z) elt ics the complex wave, where (4) 
A(z) = A(®) ew) is the spatial portion, (5) 
A(t) =a f(r) is the spatial amplitude , (6) 


Y (Zz) is the spatial phase function, 
a is the reference amplitude, 


f(F) is the spatial variation of amplitude, and 
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¥ isa position vector. 
Since photographic film measures only the time averaged intensity 


over relatively many periods of the wave, one need not be concerned with 


the temporal portion got 


of the wave. Denoting the object and refer- 
ence light beams by the subscripts o and r, ane can represent the com- 
bination of both beams at the photographic plate as: 

A= Ay+ A= Ay aa ae abe (7) 
The intensity recorded will be: 

T= Kit = (Kya K,) (KpsK) = [ad + APs AER, + ALAS) 
where * denotes complex conjugate. Assuming a linear ratio of intensity 
to transmittance of the developed hologram, the transmittance is: 

T= 1-T,2 =4-T, (lal [AI - T,(ACK.) -T.(Ae Ro) 
where T, is the transmittance corresponding to unit intensity. Upon 
reillumination by the reference beam, the transmitted wave is: 

AcT = t- Teele Al)|K, - [Tle Ko- (THAD RS (10) 
where the quantities in square brackets represent attenuation. The first 
term represents the transmission of the reilluminating reference beam, 
the second term represents the production of an attenuated object beam, 
and the third term represents a modified reconstruction of the complex 
conjugate of the object beam (Figure 3). The modification of the conjugate 
object wave depends upon the re direction of the reference beam, 


but for simple reference waveforms there is usually only a magnification 
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and angular displacement to the opposite side of the reference wave. The 
negative sign on the reconstructed wave amplitude is immaterial, since a 
negative amplitude on a continuous wave simply represents a 180° phase 
change. 

If the reillumination beam is different on the original reference 
beam, the reconstructed diffraction beams are distorted accordingly. 
For the most part, this distortion is due to different source location or 
different wavelength of light, which create magnification and angular 
displacement of the object waves. 

For laboratory use, the reconstructed holographic image can be 
considered an almost exact duplication of the original. 

3. Holographic Interferometry 

The application of holography to interferometry [Heflinger, et.al. , 
1966] derives from the fact that the holographic film records the diffrac- 
tion pattern almost linearly. When double exposed, the hologram records 
two diffraction patterns superimposed upon one another. Therefore, 
when a double exposed hologram is reilluminated, each recorded diffraction 
pattern will diffract its own first order beams. An observer will see both 
recorded scenes simultaneously. One may simply replace A, by A ol * As) 
in equations (8-10). To produce a holographic interferogram of a test 
section as in Figure 2, one exposes the hologram for one-half the exposure 


time with no flow in the test section (the uniform field provides the 
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comparison beam); then completes the exposure with the subject flow field 
present. The two reconstructed waves will interfere with each other in 
much the same manner as do the two waves of the Mach-Zehnder in the 
infinite fringe configuration. 

If, instead of making a second exposure of the hologram, one re- 
places it exactly in position and reilluminates both the noloerer and the 
object, one can again compare the two waves, except that in this case 
one of the object scenes is real and the other is a reconstructed virtual 
image. Such real-time holographic interferometry is called the "live" 
fringe technique. 

The chief advantage of either method is that, except for differences 
in the test section, both the test and the comparison beams have traveled 
through exactly the same optical regions. Optical components are then 
automatically matched. In practice rather crude optical components cant 
be used with excellent results. 

By a slight rotation of the hologram or of one of the mirrors be- 
tween exposures, one can achieve the same finite fringe presentation as 
with the Mach-Zehnder. 

a. Dark Field Interferograms 
Dark field interferograms are produced when the test section 
does not contain any light scattering properties. The observer looking at 


the reconstructed image is unable to focus on anything other than the 
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source; objects in the test section, as well as regions of interference, 
can only be observed as shadows. Such interferograms correspond to 
those taken by ordinary methods with only one direction of light beam 
passage through the test field, corresponding to only one light source. 

b. Light Field Interferograms 

True three-dimensional interferograms are obtained when a 
diffuser plate is placed between the source and the test section. The 
object beam then becomes diffused through the test section and appears 
to the observer as a continuous background of source points against which 
the test section, and its interference sea become a silouette. Be- 
cause each point of the diffuse plate acts as an individual source for any 
line of sight passing through it, equation (3) may be evaluated for any 
line which passes through both the diffuser plate and the hologram. A 
continuous evaluation of the function g can thus be provided as a function 
of position and angle. This function g, the fringe pattern, changes as 
the observer changes viewing aspect, but generally the fringe patterns 
cannot be localized by their parallax. 

By arranging several holographic plates about a test section, as 
in Figure 4, one can obtain the fringe number function for a rather con- 
tinuous segment of angular variation about the z axis of the field. The 
resulting array of integral values can be applied to the inversion of 


equation (3) to provide a solution for the density field in the test section. 
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B, THE INTEGRAL INVERSION 


1. The Basic Equation to Invert 


Equation (3), rewritten for a plane of constant z is: 


Ka 
4 (y',§ ze) =Q ys (xy 2.) dy! (11) 
where: C (xyz) = C(XY%) 24 112) 
C0 
and: Q = = (13) 
S 


x' and y' are measured in a coordinate system which is rotated by an angle 
§ about the z axis (Figure 5). 
2. Maldonado's Inversion 

The integral inversion method utilized in this investigation was 
first reported by C. D. Maldonado et.al. in 1965 [1966; Olsen, 1968]. 
It was used for obtaining plasma emissivity within a particular region 
from the measured values of emission intensity measured from outside 
the region. The form of the equation resulting from such emissivity 
studies is identical to that of equation (11). The procedure involves the 
representation of the function f of equation (11) in a complete set of 
orthogonal functions, with the expansion coefficients evaluated by use 
of the orthogonality condition. 

An earlier method involving series expansion of the function f 
had been reported in 1962 by S. I. Herlitz. The method applied to cases 
of asymmetric function of finite domain. Maldonado's method follows sim- 


ilar logic, but it is applicable to asymmetric functions of infinite domain. 
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The function f is assumed to be squared integrable over the infinite 
plane so that it may be expanded in a complete set of orthogonal functions. 
The selection of a suitable set of orthogonal functions is discussed in 
Appendix B. 

a. The Analytical Procedure 

The unknown function may be expanded in a special set of 


functions Ux eae which are orthogonal with respect to the weighting 


function er . ay). 


¥ (x,y) -y Ye los esril «UL go (KX KY) + aoe aU weel Xo) 


Wizo K=o _(o2@ 2 
‘ore 14) 


em 


where E..= & for m =0 , ©&m =1 for m =1,2,3... and Gorey a 
are the unknown complex coefficients of expansion. © is an arbitrary 
scale factor which may be considered the reciprocal of a non-dimensional- 


izing coefficient. 
tm 
The functions Uy 4,2, are defined: 
am] + 
2M K [ki Gist) |? time om,,, 
Unaoe (aX,0¥) = (-1) mt Cmeky td | © Ly ted’) (15) 


~ wk . : 
where: b = tay (2)-% , and ine is the associated 


Laguerre pol tee 
es Cag : t) = (m+ k) ae 2 ° 
_ i LTeayt cwesyr sal) CPaehy) a 


The ae Jee 


ea has a gauss transform: 


on Tm Be oy’ 


me We 16 
Tp (4758) = ju Unwaaltnaie ox’ = [k! (ask) eae : 
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where naa (xy’) are Hermite polynomials of order m+2k. 


The particular advantage of the set of functions U* is that they 


BAe 
are "invariant in form" to a rotation of coordinate system [Maldonado, 
1965]. That is, they remain an orthogonal set under a rotation of the 
coordinate system. Observe that in the equation for the polynomial (15), 
the angle ( occurs only in the complex exponential term. 


In terms of the expanded function f, and utilizing the trans- 


form relation of equation (16), equation (11) becomes: 


me 
2 M+2k 2 


q (ys) re) 5 Ev ik (mak)! 2 


M=0 k=o (17) 
6. 4 }u : 
un «)e - + Cal a) ¢ si Bn a)e Bus 


Equation (17) is subject to the following orthogonality condition: 


4iemS Fin§ 7 Pere 
) C e d§ iy H nazz XY) H 429) a 2 dy' = (18) 


3 he rk 

fous | @nsex)! (n 42x) | Wome. oven Scraislnned| 
aX 

Where 9 is the kroneker delta function. The solution of equation (18) 


applied to equation (17) leads to the expansion coefficients: 


an a [kbemak)! 
oF ee ee an (m+2k) | 
(19) 


fr 00 
SS \ 9 (y’, §)e or H saaee 7’) dy'd § 
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Equation (17), with the coefficients of equation (19), represents 
the function g by a Gram-Charlier series in the radial direction, and by 
a Fourier series in the azimuthal direction. 
Equation (14), with the coefficients of equation (19), becomes: 
0 60 ; se 
F (x,y) = a7 » 2 Ean (kt (met) VF e. (ohtreb yt) 
=o Kro 


(m +el< \ 
“Nv 00 


| (20) 
+ Real b'. 9 (y g) e — H olay) dy'dg° are errr (xx, «) 
or, inserting gaustson (lS): 


¥ (x,y) = (2 I : ea (“1 ks Oe sot) * (fx .0°y*) 


ee nek) | 
(21) 
1B (a) Pane eae (x) SIN Cmo)| € (trav 
where: 
Brrven()= : a (y' g ) COS (me ) Henson l&Y’) eine (22) 
-% -00 


oO 


AY 
ae C= sf q (y's) SIN (vn g) ene (xy') dy' d§ (23) 


Equations (21, (22) and (23) represent the fundamental 
inversion by Maldonado's method. Both analytical and numerical inversions 
were demonstrated in his series of three papers. 


The method applies quite directly to density fields. The zero 
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value of the function AG studied in interferometry is arbitrary and can 
always be chosen such that the function is zero outside of a given circular 
boundary. Singularities do not occur in real density fields except for the 
Spaces occupied by opaque objects. Cases with simple solid objects appear 
amenable to this method, but will not be discussed here. 
b. The Numerical Procedure 

Since experimental data are not obtained in analytical form, 
the indicated operations must be performed numerically. The contractual 
report resulting from the investigation by Olsen, et. al. [1966], included 
a computer program. However, it appeared more prudent to reprogram 
the basic equations to fit the available computer, and to more closely 
meet the needs of the particular application. 

The form of equation (24) is a truncated form of equation (21): 

> MK % 
£ (x,y) = (x. » - eo (sot) = e (ober 08) 


VT 
W=0 K=0 Mazk e (24) 


@ 
. [Bi (a) coslane) + Tay (x) smCma)] ECR 
The numerical evaluations of the coefficients from equations 
(22) and (23) are accomplished approximately. Since the function is 
known to be zero outside of the circular domain of Figure (5), the 
integration over y' need only be aan over finite limits, defined here as 


-S/2 to +S/2. The function g, which is always single-valued, is assumed 


to be constant within the discrete areas of the data plane represented 
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in Figures (6) and (7). This permits the integrals to be represented as 


a finite sum of integrals over the subregions: 


IMAX JAX Y34 


Beate (%) = 2 "2 (4,8 )¢ ‘costing Jas 6, He. @r)ay' @s 
G: in viet 
Wax THK Sisk von 
Draw @) DP abrs.t nS an(mng)ae Bee oo) 8 


The integrals are then evaluated over the sub-intervals with ea = 


3 
es foes Ss S 

cece = + ; Ty =" 3 , and eee = 4, 

The integrals over § are elementary, and the integrals over y' are evalu- 


ated by use of the derivative formula for Hermite polynomials [Erdely1, 


et. al. , 1953]. The a Bere formulae for B and D are: 


es, (x) = eT arta ‘ox : OARS < ) eee) 


\ =I 


| Hanvancer(*¥ 4) = One 5-4) (27) 
[= (mn S44) - cos (M44 7 


(28) 


WAX TMAX 


Dive ) = Fatma) 2 2 (i, &) 


{ acl 


. I Wonvares( V5 :)- Haveun (45-4) 


which are valid for m 7 0. The special case of m = 0 represents the axi- 
symmetric solution. The terms within the brackets [ ] are evaluated 


Oo 
by l'Hospital's rule: Bieta = ©) + and 
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ark 


lm SIN (rn Sturt )- Sw (-m ) - ¢ 
Wie Oo WM 7 t+ 
c. The Computer Program - Holofer 


In order to perform the inversions a FORTRAN IV computer 
program called HOLOFER, which includes a large number of variable para- 
meters and optional modes of operation,was written. The complete 
program, with input parameter descriptions and the working data, is 
included in appendix C. 

(1) The Basic Program. The data plane, called the G-array, 
consists of a discrete set of values of the function g of equation (11) 
taken over regular intervals of y' and 5 , of size IMAX by JMAX (Figure 7). 
The G-array is obtained in one of three modes. Mode one computes the 
G-array that corresponds to an optional test ee by the numerical 
integration of equation (11) at regular intervals of y' and § . By gener- | 
ating its own data from a test function, mode one provides the basic pro- 
gram testing capability. In mode two operation, HOLOFER interpolates 
a regular G-array from an irregularly entered set of values such as those 
normally obtained in holographic interferometry. Mode three operation 
reads in a regular array of input values directly. 

Next , equations (27) and (28) are evaluated for each set of 


coefficients B and D up to the series truncation valuesof m= MLIMIT 
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andk= KLIMIT. These B and D values are stored on external (DISK) 
storage. 

Having computed the B and D coefficients, HOLOFER next 
evaluates equation (24) for each point in the test section individually. 
Either a cartesian or a cylindrically aligned array of spatial coordinates 
can be accommodated by choice of input parameters. In the evaluation of 
equation (24), truncation of each series is performed by first specifying 
an input parameter epsilon. If a certain number of terms in a row remain 
less than epsilon times the existing partial sum, the series is truncated. 
If the index nates reach the limiting available set of B and D coefficients 
before the series has converged, the most recent few values of the partial 
= are averaged, 

(2) Add-On Functions. A second G-array may be computed and 
added to the original G-array prior to the inversion. The array is com- 
puted in the mode one operation from any given function. In addition, a 
gaussian random number of specified standard deviation may be added to 
the second G-array. The first feature provides the capability of adding 
a smoothing function to data with severe properties , such as from shock 
waves or opaque discontinuities in the interferograms. Peripheral infor- 
mation about the function may be used to determine a desirable function. 
The known add-on function is then automatically subtracted from the 


solution. The second feature allows one first to perform the inversion 
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of a set of data with a known error size and then to compare the result 
with a previous solution in order to estimate the possible error in the 
solution. 

(3) Shock Waves. The polynomial series representation of 
discontinuous or high gradient regions such as those from shock waves 
requires a large number of terms. Accordingly, the inversion routine 
requires a large index of m andk to represent steep gradients. If 
sufficiently many data points are entered on input to provide adequate 
significance in the high order terms, one can utilize rather large values 
of MLIMIT and KLIMIT to obtain reasonably accurate inversions. 

If one has prior knowledge of the function to be inverted, 
such as an isentropic solution of a flow field, an add-on function can be 
added which will smooth the net function at such discontinuities, i.e. 
the subtraction of a "plug" function shape in Figure (12). Such modifi- 
cation enhances the convergence of the resulting function, if it is well 
matched, and it therefore improves the accuracy of the result. An 
analogous method applicable to the Abel axisymmetric inversion has been 
described by Bennet, Carter, and Beredslt [1952], called the method of 
"reduced functions ,"" which operates on the fringe curve directly. 

d. Symmetry 


HOLOFER is designed to invert an asymmetric field. There 
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are, however, several simplifications to the process that result from 
symmetry in the field. 

(1) Axisymmetric Fields. For axisymmetric fields, the 
angular variations disappear with the result that only one view or one 
line of the G-array need be considered. The coefficient evaluation 
integrals over the angular domain thus become zero, for all except the 
m = 0 case, and the summation process simplifies to a single series. 
Computation time in this mode is reasonably short, from 20 - 60 seconds 
on the IBM 360-67 digital computer for KLIMIT = 1000. The majority of 
program testing has been in this mode. 

(2) With Planes of Symmetry. Planes of symmetry in the 
field reduces the angular requirements of the input data. In Figure (8), 
the data plane is shown divided into unique segments, only one of which 
need be supplied. Note that even in the asymmetric case, the data plane 
has two duplicate (though inverted) segments. The function looks the 
same from opposite sides of the domain, with reversed argument: 

2 ye.) = a (-y, & +) » Additionally, the domain of orthogonality 
is reduced for cases of planar symmetry greater than one. The regions 
which are not cross hatched, are subdomains of orthogonality. The region 
of integration for coefficient evaluation is reduced accordingly in the 
program. The index value ™M in the series representation corresponds 


to the Fourier mode of the expansion of the function f in the azimuthal 
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direction. Therefore, for cases of planar symmetry the coefficients 


are zero for certain values of m . To summarize: 


Asymmetric: uses all B and D coefficients: 

One Plane: uses B coefficients only, all 

Two Planes: uses B SALE ~.m — 0, 2,°4. 6, 
Three Planes: uses B coefficients, m=0, 3,6,9.... 
Four Planes: uses B coefficients, m=0, 4, 8,12.... 
Five Planes: uses B coefficients , m=0, 5, 10,15...... 
CEC. 


Axisymmetric corresponds to an infinite number of planes of 
symmetry: uses B coefficients m =0 only. 

While the higher numbers of symmetry planes have little practical 
value , they are included automatically in the generalization of the program 
symmetry tests. | 

3. Theoretical Results 

a. Artificial Data Generation 
Mode one of the program generates the G-array for a 
specified function. Several functions are included in subroutine FUNCT 
for use as program tests or as add-on functions. Provision is also 
made to enter an arbitrary function in subroutine SPFUN by simply 
entering the Fortran cards describing the function. In addition, a set 


of numerical values may be read in and used as a function with linearly 
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interpolated values between the given points. The G-array is evaluated 
at each value of y' and § as shown in Figure (7). Simple summation of 
the function value at ce (IMAX) points along the line y' = const.is used to 
evaluate the integral. Once the data plane is filled, the inversion process 
makes no distinction between computed or experimental data in the G- 
array. 
b. Test Inversions 

(1) ee metric. Figure (9) represents the function, com- 
puted fringe curve, and inverted solution of a gaussian function similar 
to Maldonado's test case: . =. oe x y) (30) 
The original analytical function is shown as a solid line, while the corres- 
ponding computed fringe curve is shown to a separate scale as a series of 
boxes. The greatest difference between the inverted solution and ine 
original function was 0,001, or 0.1 percent of the maximum function 
value. 

Figure (10) is a cosine-squared function: 
f= cos” (em (t+ ey 3) (31) 
The greatest inversion error was 0.5 percent of the maximum functional 
value. 
To determine the capability of the program to invert shock - 

wave type discontinuities , a circular square wave, or "plug" function,was 


inverted. Figure (11) shows a rounding off of the discontinuity that is 
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typical of a truncated series representation. The maximum error that 
occurs in the nearby region is 6.2 percent. At the center of the function, 
at radius zero, the series does not converge. The resulting 8.6 percent 
error is typical of the center point of inverted functions containing dis- 
continuities. 

Winckler , in an extensive application of the Abel inversion 
method to free jets [1948], used a hypothetical test function that is 
shown in Figure (12). Plots of the solutions he obtained are shown for 
comparison of the two methods. The characteristic "overshoot" of the 
Abel method near the discontinuity is shown, where the largest error by 
the Maldonado inversion was 3.8 percent. 

Figures (13) and (14) are inversions of typical axisymmetric 
density functions from the Winckler analysis of a free jet. They demon- 
strate the capability of the inversion routine with realistic types of | 
functions. Both inversions were accurate to within 2.6 percent. 

(2) Non-Axisymmetric. Figure (15) represents three cross- 
sections of the same gaussian function as before, now centered at 
x=.1,y=0.cm. This corresponds to the test example that Maldonado 
used in the verification of his inversion scheme. The solution is every- 
where accurate to within 0.8 percent with a very low number of series 
terms. The function is planar symmetric and thus tests only the cosine 


terms of the expansion. 
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Eioure (16) represents three-cross-sections of an asy- 
metric test case. The function is an elliptical cone of base diameters 
0.7 x 0.5cm., centered at x = 0.707, y = 0.707 cm. In the solution, the 
tip of the cone is rounded by the natural smoothing characteristics of the 
inversion method, but otherwise the maximum error is 1.5 percent. This 
test case represents a complete test of the inversion procedure, utilizing 


both B and D terms of both series expansions. 


il; EXPERIMENTAL APPROACH 


/ 


Pee DESCRIPTION OF THE APPARATUS 
1. The Holographic Table 

A holographic work table was designed and constructed to provide 
optical bench facilities about a free jet, as shown in Figures (17) and (18). 
The table was designed to evaluate various arrangements of optical com- 
ponents to achieve a wide angle of holographic field of view, and to 
suppress the vibration caused by the free jet noise. 

The table was made of pressed plywood laminate, two and one- 
half inches thick. The four by six foot table was mounted on a rotatable 
set of plywood boards, with a centerhole of eight inches diameter to 
accommodate the free jet. This arrangement provided the capability 
of rotating the entire table about the flow field. Below the rotating 


arrangement, the table rested on four small tire tubes to provide 


29 





structural vibration isolation from the building. For experiments 
conducted in C-W gas laser holography with up to thirty second exposure 
times (without the free jet) the arrangement was very successful. Be- 
neath the inner-tube mounting, the table rested on four standard auto- 
mobile type screw jacks. Recessed jack points allowed the table to be 
readily tilted to about 15° about the flow field. 

A Korad K-1 pulsed ruby laser operating at a wavelength of 6941. A. 
was employed with a Pockels cell Q switching device. The resultant 
effective hologram exposure time was about 20 ns. , eliminating problems 
with vibration during the hologram exposure. There does remain the 
problem of vibratory misalignment of the optical components between 
the two exposures of the holographic interferogram. To help damp 
accoustically-caused vibration, the mirrors were all mounted on heavy 
metal blocks. The weakest link in the setup appeared to be the beam 
splitter holders. They were lightest of the table components, and asa 
result of their vibration, holographic interferograms obtained tended to 
have a finite fringe. Occasionally a hologram would be unusable because 
the fringe spacing became too fine to resolve. 

2. The Test Section 
The table was mounted around a standing free jet. The plenum 
chamber was about 45 cm, long by 30 cm. diameter. The jet extended 


45 cm. above the plenum chamber with an inside diameter of 3.18 cm. and 
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a throat diameter of 2.0 cm. at the exit. The test section was defined 
by a square plexiglass enclosure, the four inside surfaces of which were 
inscribed with a one centimeter square grid. The grid box provided some 
vibration insulation for the jet, and also produced a self-contained co- 
ordinate system for the hologram and the corresponding photographs. 
The holograms were arranged about the grid box as shown in Figure (19), 
with the coordinate system established as shown. Commercial ground 
glass plate was used as the diffuser. 
3. Laboratory Techniques 

Alignment of the laser beam with the optical components was 
accomplished by aligning a continuous wave helium-neon gas laser through 
the rear mirror and along the ruby axis of the pulsed laser. A -20 cm. 
meniscus lens of inexpensive quality was used to diverge the beam. Care 
must be exercised in this arrangement to insure the convex side of the 
lens is towards the ruby laser. Concave surfaces in the near vicinity of 
a high power laser can themselves focus a high density of reflected energy, 
which might damage the laser itself or any near components. 

Relative intensities of reference to object ne of about 4:1 
were found to yield good holograms. 

Since the arrangement of the several mirrors, beamsplitters, 
etc. , so as to have the same optical path length from source to hologram 


for each beam is a tedious procedure , a cork pinboard was designed to 
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facilitate the table arrangements. Threads of the same length, each one 
representing a laser beam, were stretched from the laser source to their 
various holograms via different routes over a sketch of the table and held 
in position with pins. The table has a six inch grid painted on its surface 
to facilitate location of the components from the grid on the sketch. 

The pin board saved many hours of arrangement time. 

Holograms were made on Agfa-Gaefert 8E75 holographic plates. 
Development was five minutes in Kodak D-19 developer, 30 seconds in 
acetic acid stop bath of standard dilution, five minutes in Pesndees fixer, 
one minute water wash, followed by immersion in wetting agent (Kodak 
PhotoFlo) prior to drying. It was observed that good reconstructions 
of restricted field of view could be obtained immediately after developing, 
while the plate was still wet. Normal reconstructions were made with the 
continuous wave gas laser at 6328A. The resulting image magnification 
from reconstruction at a different wavelength was not considered 
deleterious. The technique used in making photographs for data extraction 


will be discussed in the next section. 


B. HOLOGRAPHIC EXPERIMENTAL RESULTS 


1. Non-Flow Holographic Experiments 


Holographic experience was gained with a Spectra-Physics model 


124 gas laser. Both back-lighted scenes and diffuse objects were recorded. 
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Real time interferometry of a test area was demonstrated and surface 
deformation interferometry was obtained. Figure (20) shows the inter- 
ference pattern on the surface of a small clamping device under strain. 
Figure (21) is one view of the interference pattern about a cigarette. Be- 
tween exposures of the hologram, the cigarette had been radiantly heating 
the platform below it for about three minutes, and the platform shows 
the interference pattern caused by the resulting thermal expansion. The 
gas laser was used to confirm the optical arrangement of Figure (19) for 
multiple hologram viewing of the free jet test scene. Successful holo- 
graphic interferograms were obtained of cigarettes at all hologram 
positions. 
2. Free Jet Experimentation 

The Korad giant pulse laser was installed on the table and holo- 
grams taken of the free jet exhausting to the atmosphere. Figure (22) 
shows a shadowgram taken at 35 psig plenum pressure. A shadowgram 
is produced directly on the hologram plate when the holographic image is 
recorded by a single exposure of the dark field technique (with the ground- 
glass diffuser absent). 

Figures (23) and (24) show two views obtained from the same 
holographic interferogram of jet flow at 60 psig. The interferometric 
data from this hologram were inverted to provide an axisymmetric 


solution. The results of the solution follow in section D. Holographic 
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interferograms were also taken at 25, 40, 45, and 50 psig. , but were not 
reduced. Turbulent flow perturbations became more predominant as the 
pressure was reduced. 

To provide a non-axisymmetric test of the method, the table was 
tilted 11° clockwise about the y axis of the table. As a result of the tilt, 
horizontal cross sections through the field became planar symmetric about 
the x-z plane. The solution of planar symmetric fields require a 90° field 
of view. Three simultaneous holographic interferograms were taken about 
the tilted jet at 60 psig. with the arrangement shown in Figure (19). Each 
of the three holograms provided a field of view of about 15°, one of which 
had several degrees obscured by the corner of the box. To provide more 
complete coverage, the table was rotated and additional holograms taken. 
Two Mentions were required. Figures (25) and (26) show the interfero- 
grams taken at 5° and 85° which were used in the data reduction. 

3. Experimental Techniques and Considerations 
Previous work at this laboratory [Sullivan, 1968] had shown the 
intensity transmission of collimated laser light beyond commercial ground 
glass falls below 30 percent of the incident intensity beyond viewing angles 
of about + 8°. Since the diffraction capability of a holographic plate ex- 
ceeds + 8°, the ground glass represents the limiting factor to the field of 


view. In fact, usable holographic interferograms were obtainable with 


from +5° to +10° field of view, centered about the object beam direction, 
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the field size being a function of the ratio of the intensities of object 
beam and reference beam. 

The holographic interferogram appears as a completely three- 
dimensional set of fringes to the observer. Only special cases of the 
fringe pattern can be localized, those corresponding to Young's fringes 
(the finite fringe background pattern), or those corresponding to regions 
of spherical or cylindrical density variations. To obtain usable data 
from the fringe number function, one must sample discrete segments 
of the information. This is done by recording a series of photographic 
interferograms at regular angles about the test section. For a completely 
general field one requires 180° field of view, while the planar symmetric 
case studied required 90°. 

There are two basic methods of obtaining sufficient angular cover- 
age of the field. The first is to take several individual holograms, 
rotating the relative angle between the hologram setup and the test 
section for each hologram. Unsteadiness in the Flow between the expo- 
sures will introduce errors. The second method involves arranging a 
series of holograms about the test section for simultaneous exposure. 
Intervals in the data from the second method are filled by interpolation. 
Interpolation over large angles requires that the function vary slowly 
in the angular direction. This experiment utilized a combination of both 


methods. 
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a. Photographic Technique 

A normal photograph of the hologram records an image of the 
focus plane as shown in Figure (28). Each position on the photograph rep- 
resents the line of sight from the image to the aperture (Figure 27). 

All of these lines of sight will represent a non-parallel set of lines. For 
reasonable camera focus distances, the deviation from parallelism is 
small and may be neglected. The spatial filtering technique shown in 
Figure (29) allows the selection of parallel sets of lines of sight for the 
recorded fringe pattern. The aperture stop at the focal plane of the lens 
filters out all but the lines parallel to the central angle. The resulting 
photographs are simpler to analyze since the angles are constant. In 
addition, fringe data from any z plane may be obtained from the single 
photograph. Mach-Zehnder interferograms, because of their single col- 
limated source, provide the same type of interferogram. 

The technique utilized for this investigation was an application 
of the lensless focusing capability of the hologram and is easier to achieve 
then the previous two methods. Figure (29a) shows a hologram being re- 
illuminated by a conjugate reference beam of small diameter. The real 
image of the test scene is formed in the same position as shown in 
Figure (3). Because the reconstruction beam is of small size, the lum- 
inated portion of the hologram represents a small aperture. The result- 


ing image has a large depth of field and a photographic film placed at the 
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image records the test scene. Additionally the rays may be focused in the 
most desirable plane by positioning the film plane, usually near the center 
of the disturbance. The lines of sight recorded represent the diverging 
bundle passing through the aperture position from the diffuser. The 
maximum angle of divergence at the edge of fe test field encountered 
was + 5°. For resolution greater than + 5°, one must compensate for the 
variation. A subroutine of the computer program was written to accom- 
plish this compensation, but the errors introduced by neglecting the bundle 
divergence have been acceptably small, and use of the routine has not been 
required. 

The specular interference of coherent light from a diffuse 
surface causes a film graininess that is inversely proportional to record- 
ing aperature size [Tanner, 1967]. Although recording apertures of one 
millimeter were used, the speckle created no problems. 

b. Data Reduction 
To obtain the photographic interferograms, a camera back 
with viewing screen was placed in the position of the real image of the 
hologram as previously described (Figure 29a). The illumination of each 
position on the hologram corresponds to a particular effective aperture 
position. The hologram was positioned such that the desired set of front 
and rear grid lines lined up on the camera back viewing screen correspond- 


ing to the desired elevation and angle of view for the picture. Reference 
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to Figures (23) and (24) show grid alignments at Z=3. cm., § =0° and 
z=1. cm., §=0° respectively. Polaroid P/N 55 film was used to record 
the image. 

To obtain the graphical plot of the fringe number function, 
the negative was used in a photo enlarger to image directly upon graph 
paper. With the enlargement adjusted to match the scale of the paper, 
the positions of the maxima and minima of the fringes along the desired 
cross section were determined visually. The fringes were counted, and 
given a graphical elevation according to their number. Initial fringe 
numbering was arbitrary, commencing with a fringe well to one side of the 
field and proceeding across the field, the light regions representing 
integer values of the fringe number. 

Figure (30) shows a typical working graph that was used for the 

§ = 5° view of the tilted jet. The finite fringe reference line was drawn 

from the field edges and the curve was transcribed. Numerical values 
were then read from the curve at regular increments of radius ( y' ) 
corresponding to the G-array size IMAX. 

For the axisymmetric case the entire procedure was repeated 
for each level z plane to be solved. Since eachaperture position on the 
hologram corresponds to a particular elevation and angle of view, a new 


picture was made for each horizontal z plane through the field. 
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For the asymmetric case, a new picture was taken for each 
angle § desired. When the angle desired was not available because of 
gaps in the holographic coverage, the curve was graphically interpolated 
' from nearby curves on either side. Fringe curves were obtained in the 
standard manner for angles on each side. The two curves were then over- 
layed on a light table and an intermediate curve drawn at the proper 
relative distance. 


IV, THE APPLICATION OF THE INTEGRAL INVERSION TECHNIQUE 
Orr eX beRIMENT AL RESULTS 


A, AXISYMMETRIC SOLUTION 

The interferometric data from the hologram of the free jet at 
60 psig. were reduced at eight positions out to z = 2.0 cm. in the manner 
described previously. Figure (31) shows the fringe curves obtained at 
“aves cimwand zZ = 1.5 cm, nae their corresponding density solutions. 
The complete set of solutions are shown topologically in Figure (32). The 
solution compares quite well with the topological features of the free jet 
solutions obtained by Ladenberg, Van Voorhis , and Winckler [1949] in their 
very extensive interferometric analysis of free jets of one centimeter 
diameter. Figure (33) is an isodensity line plot of the obtained data. 
Qualitatively, the densities compare very well with the similar plot of 


the Winckler solution (Figure 34). Exhaust density in the central region 
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of the exit plane agrees within 0.2 mg/cc (5.3 percent of the maximum 
field density) and in the central region at z/d=1. (d= nozzle diameter) 
within about 0.1 mg/cc (2.8 %), rising to a maximum of 0.4 gm/cc (10. 6%) 
between. Ata relative radius of one-half, the difference in the two 
solutions runs from about 0.05 mg/cc (1.4%) at z/d=0. rising to about 

0.2 mg/cc at z/d=0. 25 and falling to 0.1 mg/cc at z/d=1. At the jet radius, 


the maximum differences fall to about 0.05 mg/cc. 


B. NON-AXISYMMETRIC SOLUTION 

Figure (35) shows the measured fringe curves obtained for the free jet 
at 60 psig. and eleven degrees of tilt, on a horizontal plane which inter- 
sects the jet axis at a point 0.5 cm from the nozzle. A total of nine 
angular positions were sampled, every other one being shown in the figure. 
The significant trend in these data is the shoulder increasing with de- 
creasing angle on the left-hand side and with increasing angle on the right- 
hand side. A contour map of the data surface from § = 0° to § = 90° is 
shown in Figure (36). Since not all data could be obtained simultaneously , 
the run number from which the holograms used were obtained is shown. 
The non-regularities appear to be the effects of errors introduced in the 
correlation of angular views taken at different times. A data smoothing 
technique used by Maldonado, et. al. on this data plane was that of fitting 


smooth curves to the plotted points and adjusting the data to fit the 
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smoothed curves. For this experiment, however, the data were used 
without smoothing. Figure (37) shows a sketch of the tilted plane through 
the jet, with the five lines shown along which solutions were obtained. 
Figure (38) shows the solution at the five diametric cross sections from 
& = 0° to&= 90°. The effects of the shoulder variation in the input data 
show up here as variations in the position and size of the density "valley" 
which one can see in the axisymmetric topological plot of Figure (32). 
The comparison of the solution from the tilted plane at .$= 0°, 45, and 
90° is made with the solutions taken from the axisymmetric experiment 
in common plots of the two functions. The shoulder and valley features 
of the two solutions appear consistent. The central "hill," is consistent 
in the two solutions , to the point of showing a common inflection in the 
slope near radius + 0.5. The outer five points on either side of the 90° 
curve include convergence errors arising from failure of their series 
evaluations to converge. The maximum difference in the two solutions 


is about 8 percent, although the mean deviation is much less, 


C. DISCUSSION OF THE ERRORS 
1. Numerical Inversion Errors 
Several function shapes were investigated with the axisymmetric 
inversion to determine the effect of nee parameters upon the 


accuracy of the inversion. Errors resulted from two basic causes: one 
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was the failure of the series evaluation to converge within the maximum 
number of terms for which computed coefficients B and D were available; 
the other was the approximation caused by representation of the con- 
tinuous function g by a discrete set of values. 

Convergence of the series evaluation is fastest in regions where 
the functional shape matches a gaussian. As few as five terms of the 
K series were sufficient to evaluate the gaussian function to within 
one percent everywhere. The asymmetric test of Figure (15) of the dis- 
placed gaussian used only a maximum of 25 terms, corresponding to 
MLIMIT=5, KLIMIT=5 to achieve accuracy to within 0.8 percent. The 
Opposite extreme occurs when functions with steep gradients or discontin- 
uities are inverted. For example, in the test of Winckler's function 
(Figure 12), the region near the discontinuity required up to 1700 terms 
for convergence; at the discontinuity itself, convergence was not achieved 
with 2000 terms. For functions with such discontinuities, or for exper- 
imental data inversions where there are normal irregularities in the data, 
the center position (at radius zero) often fails to converge, resulting in 
error at that position of up to 10 percent. When convergence has failed, 
however , the index count of the output information clearly indicates the 
failure, and this point value may be simply neglected. Interpolation of 
the missing value from nearby convergent points normally reduces the 


error to within two percent. 
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Decreasing the number of points of the data plane reduces the 
Significance of the higher order modes of the series expansion. This rounds 
normally smooth functions, and causes oscillations in functions with dis- 
continuities. For example, the step function of Figure (11) had a maximum 
error of 6.2 percent with IMAX = 200, and maximum error of 12.7 percent 
with IMAX = 100. Additionally, the inversion at IMAX = 100 caused a severe 
oscillation over most of the step. 

Since there is no mathematical requirement in the inversion for 
regular spacing of the data grid, modification of the numerical equations 
could be made to accommodate finer mesh sizes in the regions of special 
interest. 

2. Errors in the Data 

The errors in the final solution are mostly due to errors in 
the data. There are several independent sources of data errors, covered 
here individually. 

Probably the greatest source of error in the data arises from 
the unsteadiness of the jet flow. The three apparent perturbations in 
the solution of Figure (33) correspond to 0.1 mg/cc and 0. 2 mg/cc respec- 
tively. Based upon the maximum density in the field, these correspond 
to 2.8 and 5.6 percent errors. Variations in the flow between the runs 
are considered responsible for the major fluctuations evident on the data 


plane shown in Figure (36). The asymmetric solution tends to spread 
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errors over the whole field, reducing their effect by statistical 
averaging. 

Graphical positioning of the interferogram and scale matching in 
the photo-enlarger to graph paper step is accurate to matte two percent. 
However , a two percent position error might be magnified to a five 
percent density error in regions of normally high density gradient. The 
accuracy to which one can determine the fringe number position in read- 
ing the interferogram is within + 1/8th. of a fringe, or within about 1.5 
percent of the maximum fringe number of the flow interferograms 
studied. Winckler has provided an analysis of the relative merits in 
reading several different fringe arrangements to minimize the fringe 
number error. Figure (39) shows the inverted solution of the axisymmetric 
test case at z = 0.75 cm. with a second solution superimposed. The 
second solution has been made from data with added random error of 
0.0625 standard deviation, corresponding to + 1/8th. fringe number error 
in reading the data. The resulting error varies with a maximum of 2.1 
percent. 

The background finite fringe spacing is assumed to be constant 
through the field. Figure (40) illustrates the type of error that could 
exist by the failure of linearity. An error curve that starts at zero 


at the boundaries of the flow and rises to one-half fringe in the number 
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would yield a maximum solution error of 0.18 mg/cc at the center, or 
about 5.0 percent for the axisymmetric solution. 

It is estimated in the solution presented herein that probable 
errors are less than 8 percent everywhere, corresponding Fo the maxi- 


mum difference in comparison of the two solutions in Figure (38). 


V. SUMMARY 


A. CONCLUSIONS 

A self-consistent method has been demonstrated for the acquisition 
of interferometric data about a three-dimensional density field of one 
plane of symmetry. The data set has been successfully inverted by the 
application of a recently developed mathematical inversion scheme and 
shown to be reliable to within eight percent of the density range. The 
mathematical model has been tested for realistic density patterns, 
representing supersonic wakes and jets. The pulsed laser method of holo- 
graphic interferometry can be successfully applied in environments rela- 
tively hostile to normal interferometry. There are no inherent 
restrictions to the application of the method to generally asymmetric 
fields. The optical arrangement of the system is highly flexible and can 
be modified for interferometric studies of wakes, rockets, turbo- 


machinery flows and other hostile, or highly transient events. 
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B. RECOMMENDATIONS 

Most of the errors introduced in the solution were the results of the 
data-reduction methods used. The use of an appropriately designed densito- 
meter for determining fringe locations on the negatives ne improve the 
accuracy of fringe location and merits investigation. Measurable features 
on some object in the test area should be incorporated for accurate scaling 
and positioning of the data coordinates. 

Experimental improvements could be accomplished as follows. 
An effort should be made to increase the field of view from the diffuser 
plate to provide increased angular coverage. The use of heavy, naturally 
damped mounts should be provided for the beam splitters for better 


control of the spacing of finite background fringes. 
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FIGURE 5: BASIC COORDINATE SYSTEM OF INVERSION PROCESS. 
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FIGURE 15. TEST INVERSION — DISPLACED GAUSSIAN 
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FIGURE 17. PHOTOGRAPH OF THE HOLOGRAPHIC TABLE 





FIGURE I8. PHOTOGRAPH OF THE LASER POWER SUPPLY 
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FIGURE 19. SCHEMATIC ARRANGEMENT OF HOLOGRAPHIC TABLE 
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FIGURE 24. AXISYMMETRIC FREE JET 
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FIGURE 28: EFFECT OF APPERATURE SIZE AND FOCUS PLANE 
POSITION ON PENCIL SIZE OF RAYS ABOUT A LINE- 
OF-SIGHT RECORDED BY CAMERA. 
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FIGURE 30 EXAMPLE GRAPHICAL WORK SHEET CURVES 
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FIGURE 34 SKETCH OF THE AXISYMMETRIC DENSITY SOLUTION 
OF A 1.0 CM DIAMETER FREE JET.AT 60 PSIG, 
FROM WINCKLER. 
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FIGURE 37 SKETCH OF TILTED PLANE SHOWING SOLUTION LINES. 
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FIGURE 40 ILLUSTRATING POSSIBLE ERROR DUE TO NON- 
LINEAR BACKGROUND FRINGE. 
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aE PENDIX. A 


OTHER INVERSION METHODS 


1. The Abel Inversion 
For the case of axisymmetry, the function £ becomes invariant 


to a rotation of the coordinate system, and equation eleven reduces to: 


Sf 
a(y)= 20 4 rt(r)de (Al) 
: (x2- yz)% 
where r=(x? a y?) 1/2, which is a form of Abel's integral equation with 


a well known solution. The Abel inversion has long been the standard 
method for reducing the data from axisymmetric interferograms. 
[Ladenberg, ed. 1954]. It has been well studied and optimized. Unmodi- 
fied, the Abel solution overshoots the correct values when crossing 
shock wave discontinuities. Iterative methods have been developed which 
reduce this problem [Sangster and Shaw 1968]. In addition, there is the 
"reduced function" method previously mentioned. 
2. The Generalized Abel Inversion 
Recently a generalization of the Abel inversion by A. Pal has 

been reported [Cavanaugh, et.al. 1969] which applies to the two-dimen- 
sional function of equation eleven, It is in the process of being applied 
and evaluated at another laboratory. The method appears to have very 


similar application as the method utilized in this investigation. 


79 





3. Discrete Cubes 
By representing the two-dimensional domain as a region of 
NXN, discrete cubes of constant values (see figure five), one can 
represent equation eleven as a set of N? linear equations in N? un- 
known cube densities. The complexity of the resulting matrix inversion 
problem increases with the increasing size of the Nt matrix of coeffi- 
cients , however, and the method is not presently competitive. 
4. Rowley's Method 
An analytical method involving a two-dimensional fourier trans- 
form has been reported recently [Rowley, 1969], but no application of 


the method is known to have been attempted. 
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APPENDIX B MATHEMATICAL FOUNDATIONS 


The selection of a set of functions which are orthogonal and defined 
over the domain has been made on the assumption that there exists a 
polynomial transform relationship which corresponds to the relationship 
between the functions f and g. Consider the operator notation of 
equation (11): 

Q = Sie (BI) 
which has an inverse operation: 
Fe T 9 (B2) 
One seeks a pair of orthogonal polynomial sets which, with an appropriate 
weighting function w, will satisfy the same transform operations. 
Herlitz [1963], following a procedure from the diffraction theory of 
aberrations [Born and Wolf, 1959], for the cylindrically symmetric form 
of equation (11), and for the domain r £1, used the Chebyshev functions 
[Madelung, 1953]: 
Vie(y) = SIN ( k Cosy) (B3) 


which satisfy the orthogonality relation 


+\ * 
Syn Vy), (vy) G-¥? dy S tks) sin @e)de=H Oo, (4) 


= = 


and the Zernike polynomials, which have the form of the Legendre poly- 


nomials of argument be ore | yand are orthogonal over the unit circle: 


\ 
2(2k+\) s foie (Cee = Sey (B5) 
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Using the Mehler integral [Whittaker and Watson, 1927]: 
d 


. Cos(k+a ) Yd? -  (B6) 


(Cosy - Cos ) 





MW (cosh) = 2 
T 


b YP 


where CoS ra “ and (oS po , Herlitz reports the transform 
pair: 

2A oe 7" 
and P. (4) = on Sma \/ py ) (B8) 


Thus , by expanding f(r) in the Zernike polynomials with unknown co- 
efficients 00 
Hr) =2 ARG) rs 
k=o 
and applying the transform relation (B1), one gets: 
co 


Pl 
a(y) = L — Nay) (B10) 


The orthogonality relationship (B4) allows the determination of the 


unknown coefficients 


A 
eee OL 
Application of the evaluated coefficients Ay, to the expansion of 
equation (B9) yields the desired solution for the density. 
Maldonado [1965], in a similar analysis for the axisymmetric case 
over the entire domain r< © , using a generating procedure attributed 
to Bhatia and Wolf [1954], obtained the orthogonal polynomials: 


ey (ax,ay) = (-1)° Beco cecia) (B12) 


rar 


in which Ly, are the Laguerre polynomials, defined over the interval r¢eo , 
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which have a transform relationship with gaussian weighting function: 
Sy 22 


-X X : 
: Pees 
Un. (oxy)e dx = ee Wax (ay) (B13) 
where [| ate! (x) are the Hermite polynomials. 
Hermite polynomials have the following orthogonality seecimnat: 


s Hay (ay) Hog (xy) ott y, ~ OQ ut | ade roe ak Vk (oa)] (B14) 


Expanding the function f in terms - the polynomials U: 
Zo 


2 <e oy 

E(x) =) CopCer) Ung (ox,arr) OEY) (Bis 

k-o 

and applying the transform g =T | one obtains: 
O 
7 eye (B16) 

a (y) 3 sera Cala )H, fay) e 

veo. + 


Again, the application of the orthogonality relationship (B14) provides for 


the evaluation of the unknown expansion coefficients: 
00 


Carla arene , 0 (7) Halay) dy (B17) 


As one can see, the logical structure of the development is the same as 
the Herlitz expansion, with the exception that the domain is here defined 
over the entire x, y plane. The evaluation of the unknown density function 
is made by evaluating equation (B15) with the coefficients obtained from 
B17 \, 

The same procedure is used for the asymmetric case. The poly- 


nomials chosen are [Maldonado, on 


Upaae (@xar)= (1) ‘tal (eee MEM Cts otyt) (B18) 


Wik 
which have been chosen to have the property of "invariance in form" toa 


rotation of the coordinate system [Maldonado, 1965]. 
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The "invariance in form" is demonstrated: suppose one makes a 


rotation of the axis by angle 5 . The polynomial tee becomes: 


U nor (a7) = eS rat fewny am Cina ail a Co Xx" 047) (B19) 
Or: ae 
Usnvee Kx! ay’) = e . OF. aC ay) (B20) 


Clearly, the form of the polynomial has remained unchanged. 


The polynomials , US (ax AT ) have a gaussian transform of: 
A) 
£Wy\ Re 
a 42k (§,dy') = \ Darna 2ed x 001 ) e dx 
0 


UM § Ene x* 


im, Urn Caxiare Pe" * ax’ 


+-cM E 
eC er aan (xy) 


Li Cmek 4/S ne (B21) 


where the orthogonality relationship becomes: 
A 


ty : oa 
uM en & a 
) e e 6) Eee uy) yen lay’) e a Ae 


Avi 


eee nae ere 
on Cs (me QV e \ ee Donan Oomna2yn +2.Q)\(B22) 


The property of "invariance in form" preserves the form of the relation- 
ships for both the transform operation (B21) and the orthogonality condition 


(B22) during arbitrary rotations of the axes, § . 
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A similar expansion of the function f in terms of the polynomial 


Usa iS: 
i @) o 2 
+A +N (ge ak ¥) 
ie )=) ve, Camo) Wine (nyay) e@ (B23) 


W=0 kro 
followed by a similar application of the transform operator g = Tf yields: 


Oo 0 me 
| a(m+2tg] = am 4im§ Ras 
ay, §) =) ) eq) ek mse) 2 Ce coe LH aaY’) sf 
N-0 kro (B 24) 
from which an application of the orthogonality relationship yields an 


My 
evaluation of the unknown coefficients CS ey 


( pies e 
aM gy Cone . \ a(y's) Hofer’) ume dude 


Cina die = 9/2 (mek ( 
fy 40 (B25) 


which allow the evaluation of the expanded function f of equation (B23). 
The progressive application of the same logical sequence, which has 
been given, demonstrates the development of this method of inversion 


by orthogonal polynomial expansion. 
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Pee een. 3 C 


HOROEER. 

A FORTRAN IV COMPUTER PROGRAM FOR INVERTING THE 
me tNTEGRAL DATA OF HOLOGRAPHIC INTERFEROMETRY 
POR AcvVMMETRIC FIELDS 
1. The Description of Input Parameters 

2. A Sample Terminal Dialogue 
3. The Terminal Executive Files 


4, The Program 


5. The Input Data Files Used in the Experiment 
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1. The Description of Input Parameters: 


CMS 


IMAX 


JMAX/2 


NOTE: 


KLIMIT 


MLIMIT 


KEATS 


selects CP/CMS (terminal) or OS/MVT (batch) computer 
Operation of the program. The selection is made by 
choice of the proper CMS definition card in the calling 
routine HOLOVERT. 


is the number of intervals of the input data in the radial 
dizection, grom y'-—S/2 to 45/2, 1 IMAX 201. 


1s the number of intervals of the input data in the 
azimuthal direction, from $ = 0° to the edge of the 
orthogonal subinterval. 0<JMAX/2 < 36, if JMAX/2 
is set to zero, the value of JMAX will default to 1. 
This applies to the axisymmetric case where only one 
line of data is needed. 


The dimension statements for G and GA in HOLOVERT 
set a limiting value of IMAX*JMAX at 5151. This would 
be exceeded if both values were individually set to their 
maximum. Only the dimension statements in HOLOVERT 
need be changed to change program dimensions. 


sets the number of B and D coefficients which will be 
computed corresponding to each k index of the double 
series. 

LS KLIMIT< 1000 for asymmetric (SYM = 0) 

1= KLIMIT< 2000 for symmetric (SYM _ 0) 


sets the number of KLIMIT sets of B and D coefficients 
which will be computed corresponding to each m index 

of the double series. 

1<MLIMIT. MLIMIT may be any size, limited only by 
the disk space available for storage. Values greater 
than 50 are not desirable , however. 


is the number of extra terms that will be included in 

the evaluation of the k series after the term size has 
failed the epsilon criterion. 

0 <KEXTRA€ KLIMIT. This value insures that the term 
size has stabilized at a small value prior to truncating 
the evaluation. When the evaluation fails to converge 
prior to reaching KLIMIT terms, the last KEXTRA terms 
are averaged, 
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MEXTRA the number of extra M series that will be evaluated 
after the K series evaluations have become small. 
MEXTRA has the same effect upon the m term 
evaluation as does KEXTRA on the K series. 
0<MEXTRAS MLIMIT. 


IUJEt abe the scale factor for the Hermite and Laguerre 
polynomials. 1.0<ALPHAS 5.0 
Maldonado, et.al. have shown that the factor Alpha 
has a significant effect upon the convergence 
criteria of the series. The guide line given is to 
match the edge slopes of the G function to a gaussian 
curve: Ge pee 
Then one uses the value & which matches. For the 
G functions used in this investigation, the variable 
SIZE has been adjusted to keep the function edges 
near .8 on the unit circle. ALPHA = 2.0 has been 
satisfactory for the experimental inversions. 


SIZE defines the diameter of the region from which data 
is taken, i.e. the diameter of the inversion circle. 
SIZE has been 2.5 or 3.0 for all experimental inversions , 
2.0 for all mode 1 test runs. 


EES is the convergence criterion. When the Nth term of the 
series becomes less than EPS times the partial sum to 
date, the KEXTRA or MEXTRA terms begin counting 
extra terms. If the (NtKEXTRA) term remains less 
than EPS*SUM, the series 1s truncated. EPS=.0001 
was used in most inversions. 


RHO-INF is the atmospheric density, in mg/cc. , entered as an 
input for the density conversion of the output data. 
1.176 mg/cc. was used for the experimental inversions. 


LAMBDA is the wavelength of the hologram source light in angstroms , 
\=6941.0 A for ruby laser holograms. X =6328.0 A for 
helium-neon holograins. 


BETA is the first term binomial expansion coefficient of the 
relation between index of refraction and density. 
BETA=. 000291 for red light (6900A ) BETA=. 000293 for 
yellow-green light (5500A ). 
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MODE 


SYM 


POINTS 


LINES 


DIAGNOS 


on Dewey 


PHIZER® 


selects the manner in which the G-array will be obtained. 
The three modes are described in the text and in the 
subroutine comment cards of the program. In addition 
to modes 1, 2, and 3, the parameter may be modified. 
Mode = -1, -2, or -3 will utilize subroutine FIELD2 which 
generates the B and D coefficients individually as needed. 
It requires much more time, but does not use offline disk 
space. MODE =11, 12, or 13 will assume the B and D 
coefficients are already stored on the disk, will bypass 
subroutine BDGEN and use the previously computed co- 
efficients. This provides a restart capability to the 
program. Once the coefficients have been generated, 

one can vary the output parameters for repeat inversions. 


selects the symmetry paramenters of the field to be 
inverted. SYM is called JSYM in the program. 
SYM =0 - completely asymmetric field 
Sita 125 3,4....— 1,2,3,4... regularly spaced 
planes of symmetry exist in the field, aligned with 
plane number one on the X-Z plane. 
SYM> JMAX - an axisymmetric field. 


is the number of points to be inverted on each sampling 
line. POINTS 21 


is the number of lines of points through the field to be 
inverted LINES2?1 


provides a series of levels of diagnostic printout. It is 
helpful if errors show up in the program to printout the 
parameters as they are computed. No explanation will be 
made here, as the use of this parameter involves detailed 
study of the program anyway. DIAGNOS =0. for normal 
operation. 


provides the ability to add a random error of given standard 
deviation to the add-on function G-array. The resulting 
inversion will reflect the effect of random errors to the 
data. This parameter allows some degree of error estimate 
to the experimental data. 

STD. DEV. = 0.0 for normal operation 

STD. DEV. = .0625 for + 1/8 fringe number error in the data. 


is the angle in degrees, which the first line to be inverted 
will make with the X axis. PHIZERO =0.0 for most uses. 
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DELPHI is the angular incrementation, in degrees, that successive 
inversion lines will make with the X-axis. The second line 
will be at PHI=PHIZERO+DELPHI, etc. 


YPZERO is the y' position of the first line to be inverted, in the 
primed coordinate system rotated by angle PHI 


YPRANGE js the range over which y' will uniformly vary for the 
number of lines to be inverted. If one sets DELPHI=0. 
and YPRANGE #40. A cartesian array of points will be 
selected, if DELPHI#0. and YPRANGE=0. A cylindrical 
array of points will be selected. 


XPZERO is the position of the beginning of each line to be inverted 
(in the primed coordinate system). 


XPRANGE is the range over which the number of points selected will 
range, beginning at XPZERO. 


TST.FUN is the number of the function to be inverted. In mode 1 
this number corresponds to the function number in sub- 
routine FUNCT. In modes 2 and 3 it serves as a run 
identifying number on the output. 


ADD.FUN is the number of the function in FUNCT to be used as an 
add-on function. ADD. FUN=0. will default the use of an 


add-on function. 


GARRAY is the output parameter for printout of the G-array used. 


GARRAY = 
0, -- no printout. 
1. -- prints the G-array. 
2. -- prints the G-array with a pause for paper spacing 


on the terminal prior to printing. 
3. -- prints both G-array and the add-on G-array. 


4, -- prints both G-arrays with terminal pauses for 
spacing. 

-1, -- punches the G-array on a deck of cards and 
prints the G-array 

-3, -- punches the modified G-array sum of both on 
a deck of cards, then prints both G-arrays 
separately. 
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GRAPH 


LIN. PRT 


PAP, BND 


prints a graph of each line of the G-array. GRAPH = 


Ono Crap. 

1. -- a graph of G-array for each J. 

2. -- same as 1., with terminal pause. 

3. -- both G-array and add-on G-array are on each graph. 
4, -- same as 3. , with terminal pause. 


prints the inversion results LIN. PRT = 


Q. -- G-array is set up, inversion results are not printed. 
Inversion is not performed unless MAP. BND=0. 


1, -- prints the inversion results a line at a time. 

2. -- prints the inversion results with a terminal pause 
prior to each line of points, for paper spacing. 
It also sets a terminal pause prior to MAP output, 
if MAP is used. 


prints a map of the inversion array on the printer. Uses 
library subroutine MTMPII. MAP.BND= 


0. -- map is not printed. 
#0. -- the number used defines the interval between con- 
tours of the function f on the map, e. g. 0.2. 


A,B,C,D,E,P are variable function parameters for the functions in 


Subroutine FUNCT, as used in mode I inversion. 


S,T,U,V,W,Q are variable function parameters for the functions of sub- 


PHisyM 


KO; YO 


routine FUNCT, as used for the add-on function. 


is an input angle entered in the data cards for mode 2 or 3 
operation that defines a rotation between the inversion 
coordinate system and the laboratory, or output coordinate 
system, e.g. if the plane of symmetry of the function was 

at 5° in the laboratory system, one would enter 5.0 in the 
data cards. Output would then be in the laboratory coordinate 
system. Raw data for mode 2 is entered in the laboratory 
coordinate system; G-array data for mode 3 is entered in the 
data plane coordinate system. 


are the laboratory coordinates of the estimated center of the 
function, entered in the data cards used in mode 2, to provide 
a well-centered inversion. Output is in laboratory coordinates. 
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INPUT DATA for batch operation, input data cards are read in order: 


1) 


7 cards of input data, format 89 of calling program 
HOLOVERT; data corresponds exactly to parameters 
in the sample dialogue. 


IF USED: data cards for numerical function, formats 
88, 89 of subroutine FREAD. 


IF USED: data cards for numerical add-on function, 
formats 88, 89 of subroutine FREAD. 


input data for mode 2 or mode 3 


-mode 2: formats 59, 58 of subroutine SHEET; 
if NCODE.GE.1: format 29 of subroutine 
SIM. 

-mode 3: formats 39, 38 of subroutine READ. 


a blank card if another data set follows, otherwise a 
FINISH card; format 60, HOLOVERT. 


For terminal operation, input data 2), 3), and 4) are 
entered in the above formats in separate files and 
given optional filenames/filetypes. The organization 


of files for execution is handled by exec file HOLOFER. 


A sample terminal dialogue follows. 
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